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We present precise predictions for Higgs boson production in association with a jet. Our calcula¬ 
tion is accurate to next-to-next-to-leading order (NNLO) QCD in the Higgs Effective Field Theory 
and constitutes the first complete NNLO computation for Higgs production with a final-state jet in 
hadronic collisions. We include all relevant phenomenological channels and present fully-differential 
results as well as total cross sections for the LHC. Our NNLO predictions reduce the unphysical 
scale dependence by more than a factor of two and enhance the total rate by about twenty per¬ 
cent compared to NLO QCD predictions. Our results demonstrate for the first time satisfactory 
convergence of the perturbative series. 


Further exploration of the Higgs boson discovered by 
the ATLAS and CMS collaborations [1, 2] will be a pri¬ 
mary goal of the continued experimental program of the 
LHC. In the well-measured decay modes, h —► 77, WW 
and ZZ , the determined couplings agree with the Stan¬ 
dard Model (SM) values to 20 — 30 percent [3, 4]. This 
agreement will be further probed during the upcoming 
LHC run. The predictions of the SM are expected to be 
tested to the five percent level in several production and 
decay modes [5] , providing an exciting opportunity to dis¬ 
criminate between different mechanisms of electroweak 
symmetry breaking. In addition, new properties of the 
Higgs boson will be accessed through the measurement of 
its kinematic distributions. These measurements will test 
whether the tensor structures of the Higgs couplings are 
correctly predicted by the SM, whether additional oper¬ 
ators mediate Higgs production and decay, and whether 
there are new particles that modify the loop-induced ggH 
and 77# couplings. 

Such studies [6] are currently limited by the available 
statistics. However, this situation will change during Run 
II of the LHC, and eventually the limiting factor in the 
search for deviations in Higgs boson properties from pre¬ 
dictions will be our understanding of SM theory. This is 
apparent from a recent coupling combination performed 
by ATLAS [7]. The uncertainty on the theoretical predic¬ 
tions dominates the systematic error budget in all of the 
di-boson decay modes. One component of this error is 
the overall signal normalization, for which a precise cal¬ 
culation of inclusive Higgs production in the gluon-fusion 
channel is needed. The completion of the next-to-next- 
to-next-to-leading order (N 3 LO) calculation was recently 
announced [8]. The other major component of the theo¬ 
retical error is the need for improved predictions of the 
differential spectra that enter every analysis. In some 
final states this need is obvious; for example, in the di- 
leptonic decay of the WW channel a mass peak cannot 
be reconstructed, and accurate calculations of both signal 
and background distributions are needed in order to de¬ 
vise appropriate experimental search strategies. Even in 


modes where a resonance peak can be reconstructed, such 
as the 77 channel, the Higgs candidate events are cate¬ 
gorized according to their transverse momentum (p±) in 
order to improve the signal significance. Such a division 
relies upon accurate and precise theoretical predictions 
for the Higgs p± and rapidity distributions, that are used 
to reweight the parton-shower Monte Carlo simulations 
employed by the experimental collaborations. 

In this Letter we take a major step toward improv¬ 
ing SM predictions for several kinematic distributions 
employed in the analysis of Higgs boson properties, by 
providing a next-to-next-to-leading order (NNLO) calcu¬ 
lation of Higgs boson production in association with a 
jet. Compared to previous computations [9, 10], we in¬ 
clude all relevant channels and color structures. The phe¬ 
nomenological impact of this result spans all Higgs search 
channels. In the WW final state it refines the division of 
the signal prediction into exclusive zero-jet, one-jet and 
inclusive two-jet bins, and it can be used to improve the 
resummation of the jet-veto logarithms that accompany 
this division [11], For all final states our calculation can 
be used to more accurately re-weight the Higgs p± distri¬ 
bution obtained from Monte Carlo. Finally, it will allow 
for the comparison of the measured differential distribu¬ 
tions from LHC Run II with the most precise SM theory 
to more incisively probe the mechanism of electroweak 
symmetry breaking. 

Our calculation also represents a technical milestone 
in the application of perturbative QCD to the modelling 
of hadronic collisions. The past few years have seen a re¬ 
naissance in the development of subtraction techniques 
designed to turn our knowledge of the infrared struc¬ 
ture of QCD at NNLO into actual phenomenological pre¬ 
dictions for hadron-collider observables [9, 10, 12]. Our 
result demonstrates the power of these newly-developed 
methods in assisting the continued exploration of Nature 
at the LHC. 

Our Letter is organized as follows. We first review the 
theoretical framework that we use to obtain the results 
reported in this paper. Since this has been described in 
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Figure 1: Cancellation of 1/e poles in the qg channel. Note Figure 2: Dependence of the total LO, LO and NNLO cross- 

that individual contributions have been rescaled by a factor sections on the unphysical scale fi. See text for details, 

of 0.1, while the sum of them is not rescaled. 


detail in our previous work on Higgs plus jet production 
in pure gluodynanrics [9] , we only sketch here the salient 
features of the calculation. We then present the numer¬ 
ical results of the computation including NNLO results 
for cross sections of Higgs plus jet production at various 
collider energies and for various values of the transverse 
momentum cut on the jet. We also discuss the NNLO 
QCD corrections to the transverse momentum distribu¬ 
tion of the Higgs boson. Finally, we present our conclu¬ 
sions. 

We begin by reviewing the details of the computation. 
Our calculation is based on the effective theory obtained 
by integrating out the top quark. For values of the Higgs 
p l below 150 GeV, this approximation is known to work 
to 3% or better at NLO [13, 14]. Since the Higgs boson re¬ 
ceives its transverse momentum by recoiling against jets, 
we expect that a similar accuracy of the large-mt ap¬ 
proximation can be expected for observables where jet 
transverse momenta do not exceed 0(150) GeV as well. 

The effective Lagrangian is given by 

c = ~\g^G^ v +£ qopqi - ( 1 ) 

i 

where Gj?J is the gluon field-strength tensor, H is the 
Higgs boson field and q r denotes the light quark field 
of flavor i. The flavor index runs over the values i = 
u, d , s, c, 6, which are all taken to be massless. The co¬ 
variant derivative dp contains the quark-gluon coupling. 
The Higgs vacuum expectation value is denoted by v, 
and Ci is the Wilson coefficient obtained by integrating 
out the top quark. The calculation presented here re¬ 
quires C\ through O(a^), which can be obtained from 
Ref. [15]. Both the Wilson coefficient and the strong 
coupling constant require ultraviolet renormalization; the 
corresponding renormalization constants can be found 
e.g. in Ref. [16]. 

Partonic cross sections computed according to the 
above prescription are still not finite physical quantities. 


Two remaining issues must be addressed. First, contribu¬ 
tions of final states with different number of partons must 
be combined in an appropriate way to produce infrared- 
safe observables. This requires a definition of final states 
with jets. We use the anti -hr jet algorithm [17] to com¬ 
bine partons into jets. Second, initial-state collinear sin¬ 
gularities must be absorbed into the parton distribution 
functions (PDFs) by means of standard MS PDF renor¬ 
malization. A detailed discussion of this procedure can 
be found in Ref. [18]. 

The finite cross sections for each of the partonic chan¬ 
nels ij obtained in this way have an expansion in the MS 
strong coupling constant a s = a s (/i), defined in a theory 
with five active flavors, 
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Here, the omitted terms indicated by 0{a\) include the 
af factor that is contained in the leading order cross sec¬ 
tion crj°). Our computation will include the gg and qg 

partonic cross sections at NNLO, a gg and a qg . where q 
denotes any light quark or anti-quark. At NLO, it can be 
checked using MCFM [19] that these channels contribute 
over 99% of the cross section for typical jet transverse 
momentum cuts, p± ~ 30 GeV. We therefore include the 
partonic channels with two quarks or anti-quarks in the 
initial state only through NLO. 

In addition to the ultraviolet and collinear renormal¬ 
izations described above, we need the following ingre- 
(2) (2) 

dients to determine a gg and a q " g : the two-loop vir¬ 
tual corrections to the partonic channels gg —4 Hg and 
qg —4 Hq; the one-loop virtual corrections to gg —> Hgg , 
gg —> Hqq and qg —>• dlqg ; the double real emission 
processes gg -> Hggg , gg Hgqq , qg -> Hqgg and 

qg —y HqQQ 1 where the QQ pair in the last process can 
be of any flavor. The helicity amplitudes for all of these 
processes are available in the literature. The two-loop 
amplitudes were computed in Ref. [20] . The one-loop cor¬ 
rections to the four-parton processes are known [21] and 
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available as a Fortran code in the MCFM program [19]. 
For five-parton tree-level amplitudes, we use compact re¬ 
sults obtained using BCFW recursions [22]. 

The difficulty in completing the NNLO calculation be¬ 
comes apparent when one attempts to combine these con¬ 
tributions and cancel the infrared divergences that ap¬ 
pear separately in each component. The problem is that 
final states with different multiplicities live in different 
phase-spaces; this feature makes it impossible to combine 
them directly. The issue becomes obvious if one looks at 
how 1 /e singularities appear in different contributions. 
Indeed, the 1/e poles coming from loop amplitudes are 
explicit ones, but those coming from the real-emission 
corrections only appear upon integration over the unre¬ 
solved region of phase space. However, since we want to 
keep the calculation fully differential, we want to avoid 
integrating over the phase-space for higher-multiplicity 
processes. 

To reconcile these two requirements, which at first 
sight appear to be mutually exclusive, we use the sector- 
improved residue subtraction approach [23-26]. This is 
an outgrowth of the sector-decomposition method [27-29] 
used to compute the differential cross sections for Higgs 
boson and electroweak gauge bosons through NNLO. Sec¬ 
tor decomposition uses the observation that the relevant 
singularities can be isolated using appropriate pararne- 
terizations of phase space and expansions in plus distri¬ 
butions. Sector-improved residue subtraction combines 
this with the idea that a pre-partitioning of the final- 
state phase space similar to the FKS subtraction used at 
NLO [30] allows to extend this technique to 2 —> 2 and 
more complicated scattering processes. A detailed dis¬ 
cussion of the phase-space parameterizations needed to 
handle all the contributing partonic processes was given 
in Ref. [9], to which we refer the reader for more details. 
Note, however, that Ref. [9] dealt with gg —> H + ng par- 
tonic processes for which both the phase-space and the 
matrix elements are highly symmetric. For the quark- 
gluon channel this symmetry is lost and one has to con¬ 
sider a larger number of “sectors” compared to the case 
of pure gluodynamics. 

Before discussing numerical results, we would like to 
point out two things in connection with the application 
of sector-improved residue subtraction method. First, 
we note that upon applying this method, one automat¬ 
ically generates subtraction terms that allow extraction 
of 1 /e singularities and, at the same time, make integra¬ 
tion of the finite remainders possible. The key point is 
that these subtraction terms are obtained from univer¬ 
sal soft and collinear limits of scattering amplitudes that 
were computed long ago in Refs. [31-38]. The univer¬ 
sality of these subtraction terms makes the method of 
improved-sector decomposition attractive and, in princi¬ 
ple, applicable to processes of arbitrarily high multiplic¬ 
ity. Second, when sector-improved residue subtraction is 
applied to a physical process, it leads to a Laurent ex¬ 


pansion of the various contributions to the cross section 
in the dimensional regularization parameter e; the coeffi¬ 
cients of this expansion are computed numerically. Since 
final physical cross sections are independent of the reg¬ 
ularization parameter, the quality of l/e n , n = 4, 3,2,1, 
cancellation is an important check of the correctness of 
the implementation of the method. To show the quality 
of the cancellation in our case, in Fig. 1 we present var¬ 
ious contributions to the 1 /e pole of the partonic cross 
section, together with the residual non-cancellation, in 
the qg channel. We show these quantities as functions 
of the distance from the partonic threshold, defined as 
P = y/1 - Sth/s, y/s^h = yjmjj +p 2 j_ yCUt +p±,cut ■ We 
see that the cancellation is very good, at the level of one 
per mill or better. Although in Fig. 1 we display the re¬ 
sult for the total cross section, we have also checked that 
the cancellation holds at a similar level for kinematic dis¬ 
tributions. 

In addition, we have checked that our results for 
the Higgs plus 2-jet cross section at NLO agree with 
MCFM [39], for both the fiducial cross section and 
for several kinematic distributions. We have two sepa¬ 
rate numerical implementations of the sector-improved 
residue subtraction method that demonstrate good 
agreement. Furthermore, an independent calculation 
was also performed using the jettiness-subtraction tech¬ 
nique [40] , and good agreement for the fiducial cross sec¬ 
tions was found. 

We now turn to the discussion of numerical results. 
We first compute the LO, NLO and NNLO cross sec¬ 
tions for Higgs plus jet production pp —► H + j at the 
8 TeV and 13 TeV LHC. We use rn H = 125 GeV and 
m t = 172.5 GeV. To define the cross section, we use the 
anti-/cj_ algorithm with A R = 0.5 and a cut on the jet 
transverse momentum p± > 30 GeV. We employ parton 
distribution functions and the strong coupling constant 
as provided by the NNPDF21LO [41], NNPDF23NLO 
and NNPDF23NNLO [42] PDF sets to compute respec¬ 
tively LO, NLO and NNLO cross sections. We set the 
renormalization and factorization scales to the mass of 
the Higgs boson and we estimate the uncertainty associ¬ 
ated with higher orders in perturbation theory by chang¬ 
ing the scale by a factor of two. For the 8 TeV LHC, 
we find a pp ^H+j = 3.9/ji.i pb, 5.pb, 6.7/[/g pb at 
leading, next-to-leading and next-to-next-to-leading or¬ 
der, respectively. Results for p = m#/2 and p = 2to# 
are shown as super- and sub-scripts, respectively. For 
p = run, the NLO (NNLO) cross section exceeds the 
leading order one by 44% (72%), indicating reasonable 
convergence of perturbative expansion. The convergence 
is better for lower scales: for example, for p = itih/2 the 
NLO (NNLO) cross section exceeds the leading order one 
by 23% (29%). As expected, the scale uncertainty is sig¬ 
nificantly reduced at NNLO. This is also illustrated in 
Fig. 2, where we plot the total cross section at LO, NLO 
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Figure 3: Higgs plus jet production cross-sections in depen¬ 
dence of the cut on the jet transverse momentum. The mini¬ 
mal cut we consider is p± > 30 GeV. See text for details. 


and NNLO as a function of the unphysical scale p over 
the range p £ [ p±, C ut '■ 2 mu}- We estimate the residual 
uncertainty due to PDF to be at the 0(5%) level. The 
situation is similar for the 13 TeV LHC. More precisely, 
we find a pp ^ H+ j = 10 . 2 + 2'.6 Pb> 14.7i|;5 Pb> 17-5pb 
at leading, next-to-leading and next-to-next-to-leading 
order, corresponding to a NLO (NNLO) increase with re¬ 
spect to LO of 44% (72%) for p = run and of 25% (31%) 
for p = m#/2. 

It is interesting to understand to what extent pertur¬ 
bative QCD corrections depend on the kinematics of the 
process and/or on the details of the jet algorithm. One 
way to study this is to explore how the NNLO QCD cor¬ 
rections change as the lower cut on the jet transverse mo¬ 
mentum is varied. We show corresponding results for the 
8 TeV LHC in Fig. 3 where the cumulative distribution 
for cr(H+j,p± j > pj_,cut) is displayed. The inset in Fig. 3 
shows ratios of NNLO(NLO) to NLO(LO) H + j cross- 
sections, respectively, computed for pp = Hr = tuh as 
a function of the jet pj_-cut. It follows from Fig. 3 that 
QCD radiative corrections depend on the kinematics. In¬ 
deed, the NNLO to NLO cross-sections ratio changes 
from 1.25 at p± = 30 GeV to ~ 1 at p± ~ 150 GeV. 

In Fig. 4 we show the Higgs boson transverse momen¬ 
tum distribution in the reaction pp —> H + j, for three 
consecutive orders of perturbation theory. We require 
that there is a jet in the final state with a transverse mo¬ 
mentum higher than p± j > 30 GeV. Note that the two 
bins closest to the boundary p± t H = 30 GeV have been 
combined to avoid the well-known Sudakov-shoulder ef¬ 
fect [43]. Away from that region, the NNLO QCD radia¬ 
tive corrections increase the NLO cross-section by about 
20%, slowly decreasing as p±,h increases. 

In conclusion, we have presented a calculation of the 
NNLO QCD corrections to the production of the Higgs 
boson in association with a jet at the LHC. This is the 
first complete computation of NNLO QCD corrections to 
a Higgs production process with a jet in the final state. It 
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Figure 4: Higgs boson transverse momentum distribution in 
pp —> H+j at 8 TeV LHC. The jet is defined with the anti-fcx 
algorithm with A R = 0.5 and the cut on the jet transverse 
momentum of 30 GeV. Further details are explained in the 
text. 


shows that techniques for performing NNLO QCD com¬ 
putations, that were in the development phase for several 
years, can indeed be used to provide precise predictions 
for complex process at hadron colliders. The total cross 
section for if + jet production receives moderate NNLO 
QCD corrections. For jets defined with the anti-fcj_ algo¬ 
rithm with p±j > 30 GeV, we find NNLO QCD correc¬ 
tions of the order of 20% for p = mn■ These moderate 
corrections are the result of the smaller corrections for 
the qg channel w.r.t the gg one, and a suppression of the 
gg channel due to qq final states not considered in previ¬ 
ous analyses [9, 10]. Beyond the total cross section, our 
computation will have important implications for many 
processes that are used to study properties of the Higgs 
boson, including W + W~ and 77 final states, primarily 
through improved modelling of the Higgs transverse mo¬ 
mentum and rapidity distributions. In particular, since 
the complete N 3 LO computation of the Higgs boson pro¬ 
duction cross section is available, a consistent computa¬ 
tion of the 77-1-0 jets, H + 1 jet, H + 2 jet and H + 3 jet 
exclusive processes becomes possible for the first time. 
Furthermore, since the Higgs boson is a spin-zero parti¬ 
cle, our computation can be easily extended to include 
Higgs boson decays, to enable theoretical predictions for 
fiducial cross sections and kinematic distributions for the 
particles that are observed in detectors. Once this is 
done, our calculation will provide a powerful tool that 
will help to understand detailed properties of the Higgs 
boson at the LHC. 
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